cap log close 
log using "${logpath}CS_estimator.log", replace

/*******************************************************************************
CS_estimator.do

Analyzes district-level data using the Callaway and Sant'anna estimator
*******************************************************************************/

qui {
	use "${cleandatapath}dist_data.dta", clear
	drop if dist_code==14 // Dropping Connecticut

	foreach var in treat dropme plotnum point se ul ll {
		gen `var'=.
	}
}

local pre=24
local post=23

foreach var in PA civ crim {
	qui {
		replace dropme=0
		if "`var'"=="crim" {
			local xti = "Months Since First Criminal Enforcement"
			local outlist "log_tot_pmt log_rides"
		}
		if "`var'"=="civ" {
			local xti = "Months Since First Civil Enforcement"
			local outlist "log_tot_pmt log_rides"
		}
		if "`var'"=="PA" {
			local xti = "Months Since Prior Authorization Implementation"
			replace dropme=1 if month_date<612
			local outlist "tot_pmt log_tot_pmt ambulance_id log_rides active_firms log_active_firms"
		}
		replace treat=mofd(`var'_date)
		replace treat=0 if missing(treat)
		
		su month_date if dropme==0
		local firsttime=`r(min)'
		local lasttime=`r(max)'
		replace dropme=1 if treat!=0 & treat<`firsttime'+`pre'-1
		replace dropme=1 if treat!=0 & treat>`lasttime'+`post'
		
		su treat if treat!=0 & dropme==0
		local pre_coeffs=`r(max)'-`firsttime'
		local post_coeffs=`lasttime'-`r(min)'
		local total_coeffs=`pre_coeffs'+`post_coeffs'
	}

	forvalues cont=0/1 {
		if `cont'==0 {
			local notyet = ""
		}
		if `cont'==1 {
			local notyet = "notyet"
		}
			
		matrix `var'_res=[0,0,0,0,0,0,0,0,0]'

		if "`var'"!="PA" | `cont'==0 {
			foreach out in `outlist' {
				qui {
					csdid `out' if dropme==0, ///
						ivar(dist_code) time(month_date) gvar(treat) agg(event) long `notyet'
					estat cevent, window(0 `post')
					matrix `var'_res=[`var'_res, r(table)]
				
				*Note csdid, long returns difference between t-1 and relevant period for the pre-period, we want the negative of this:
				}
				
				if "`out'"=="log_tot_pmt" {
					qui {
						replace plotnum=_n if _n<=`total_coeffs'
						replace point=e(b)[1,plotnum]
						replace se=sqrt(e(V)[plotnum,plotnum])
						replace ul=point+se*1.96
						replace ll=point-se*1.96
						replace plotnum=plotnum+1 if _n>=`pre_coeffs'
						replace plotnum=`pre_coeffs' if _n==`total_coeffs'+1
						replace point=0 if plotnum==`pre_coeffs'
						foreach bound in point ul ll {
								replace `bound'=-`bound' if plotnum<`pre_coeffs'
						}
						replace plotnum=plotnum-`pre_coeffs'+`pre'
					}
					
					sort plotnum
					if (`cont'==0) | (`cont'==1 & "`var'"!="PA") {
						twoway (connected point plotnum) ///
							(rcap ul ll plotnum, color(navy)) ///
							if plotnum>=1 & plotnum<=`pre'+`post'+1, ///
							xti("`xti'") ///
							xline(25, lcolor(gs8) lpattern(dash)) ///
							xlab(1 "-24" 7 "-18" 13 "-12" 19 "-6" 25 "0" ///
								31 "6" 37 "12" 43 "18" 49 "24") ///
							yti("Total Payments in District (Logs)") ///
							yline(0, lcolor(gs8)) ylab(, angle(0)) ///
							leg(off) graphregion(color(white)) bgcolor(white)
						if "`var'"=="PA" {
							qui graph export ${outpath}Figure_A4a.pdf, as(pdf) replace // Figures A4a
						}
						if "`var'"=="civ" & `cont'==0 {
							qui graph export ${outpath}Figure_A4c.pdf, as(pdf) replace // Figures A4c
						}
						if "`var'"=="crim" & `cont'==0 {
							qui graph export ${outpath}Figure_A4d.pdf, as(pdf) replace // Figures A4d
						}
						if "`var'"=="civ" & `cont'==1 {
							qui graph export ${outpath}Figure_A7a.pdf, as(pdf) replace // Figures A7a
						}
						if "`var'"=="crim" & `cont'==1 {
							qui graph export ${outpath}Figure_A7b.pdf, as(pdf) replace // Figures A7b
						}
					}
					foreach repvar in plotnum point se ul ll {
						qui replace `repvar'=.
					}
				}
				if "`out'"=="log_active_firms" & `cont'==0 {
					qui {
						replace plotnum=_n if _n<=`total_coeffs'
						replace point=e(b)[1,plotnum]
						replace se=sqrt(e(V)[plotnum,plotnum])
						replace ul=point+se*1.96
						replace ll=point-se*1.96
						replace plotnum=plotnum+1 if _n>=`pre_coeffs'
						replace plotnum=`pre_coeffs' if _n==`total_coeffs'+1
						replace point=0 if plotnum==`pre_coeffs'
						foreach bound in point ul ll {
								replace `bound'=-`bound' if plotnum<`pre_coeffs'
						}
						replace plotnum=plotnum-`pre_coeffs'+`pre'
					}
					
					sort plotnum
					twoway (connected point plotnum) ///
						(rcap ul ll plotnum, color(navy)) ///
						if plotnum>=1 &	plotnum<=`pre'+`post'+1, ///
						xti("`xti'") ///
						xline(25, lcolor(gs8) lpattern(dash)) ///
						xlab(1 "-24" 7 "-18" 13 "-12" 19 "-6" 25 "0" ///
							31 "6" 37 "12" 43 "18" 49 "24") ///
						yti("Active Firms in District (Logs)") ///
						yline(0, lcolor(gs8)) ylab(, angle(0)) ///
						leg(off) graphregion(color(white)) bgcolor(white)
					qui graph export ${outpath}Figure_A4b.pdf, as(pdf) replace // Figure A4b
					foreach repvar in plotnum point se ul ll {
						qui replace `repvar'=.
					}
				}
			}
		
			if "`var'"=="PA" {
				qui putexcel set "${outpath}Table_A2.xlsx", replace
				matrix rownames `var'_res = b se z pvalue ll ul df crit eform
				matrix colnames `var'_res = blank tot_pmt log_tot_pmt ambulance_id log_rides active_firms log_active_firms 
			}
			else if `cont'==0 {
				qui putexcel set "${outpath}Table_A3_`var'.xlsx", replace
				matrix rownames `var'_res = b se z pvalue ll ul df crit eform
				matrix colnames `var'_res = blank log_tot_pmt log_rides 
			}
			else {
				qui putexcel set "${outpath}Table_A14_`var'.xlsx", replace
				matrix rownames `var'_res = b se z pvalue ll ul df crit eform
				matrix colnames `var'_res = blank log_tot_pmt log_rides 
			}
			qui putexcel A1=matrix(`var'_res), names // Tables A2, A3, and A14
		}
	}
}

log close